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Abstract 

Von Hippel-Lindau (VHL) syndrome is a hereditary condition predisposing to the development of different cancer forms, 
related to germllne Inactivatlon of the homonymous tumor suppressor pVHL. The best characterized function of pVHL is the 
ubiquitinatlon dependent degradation of Hypoxia Inducible Factor (HIF) via the proteasome. It Is also involved In several 
cellular pathways acting as a molecular hub and Interacting with more than 200 different proteins. Molecular details of 
pVHL plasticity remain in large part unknown. Here, we present a novel manually curated Petri Net (PN) model of the main 
pVHL functional pathways. The model was built using functional Information derived from the literature. It includes all major 
pVHL functions and Is able to credibly reproduce VHL syndrome at the molecular level. The reliability of the PN model also 
allowed In silico knockout experiments, driven by previous model analysis. Interestingly, PN analysis suggests that the 
variability of different VHL manifestations is correlated with the concomitant inactivatlon of different metabolic pathways. 
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Introduction 

Pathological deregulation of cellular pathways often results in a 
family of complex and correlated diseases commonly termed 
cancer [1]. Cancer is a multi factorial disease where different 
causes contribute to its development. Several computational 
methods have been developed to explore the functional pathways 
involved in tumorigenesis. Some of them focus on differential gene 
expression between healthy and pathologic tissues [2-3], on 
protein-protein interaction network analysis [4-5] or on molecular 
dynamics simulations [6]. Other methods approach the disease 
through discretization of pathological components that result in 
tumor [7]. All of these approaches are very powerful when the 
variables related to the disease, although complex, are well known 
and studied. A multi-factorial disease can be approached by means 
of mathematical theory, building a theoretical model where cell 
components are connected with each other. In biology, several 
problems were dealt with network theory [8-9]. A network is a 
group of objects strongly inter-connected with each other (e.g. 
proteins and enzymes of a pathway or animals belonging to 
interacting populations). Their construction and subsequent 
simulation is made via mathematical analysis of the connections 
between nodes found in the system and their time-dependent 
behavior [10]. A biological network is generally composed of 
proteins, nucleic acids and cofactors connected by biological 
reactions such as protein complex formation or enzyme activity 
regulation [10]. Von Hippel-Lindau syndrome (VHL) [11] is a 



good study case to test the network theory apphed to cancer due to 
the similar medical history and pathological phenotype that 
patients share. WhUe hereditary cancers represent only a small 
part of all human tumors, their investigation represents a challenge 
to understand the pathway leading to tumor formation. In 2010, 
Heiner et al. first approached VHL using the so-called Petri Net 
(PN) simulation networks [12]. Their work, inspired by a previous 
theoretical model of cellular oxygen-related pathways [13] [14], 
was a preliminary investigation of the core oxygen sensing system 
and its connection with VHL onset. Heiner and coworkers 
proposed three different functional modules responsible for 
hypoxia network control and for HIF- lot degradation [12]. In 
other words, they theorized that hereditary forms of cancer, such 
as different manifestations of VHL, are the result of different and 
concomitantly compromised metabolic pathways. 

Von Hippel-Lindau Disease 

Von Hippel-Lindau protein (pVHL) is the product of the von 
Hippel-Lindau gene, located in the short arm of 3rd chromosome, 
and constantly transcribed in both fetal and adult tissues [15]. 
Mutations of pVHL are related to a pathological outcome termed 
VHL syndrome, an inherited form of cancer [16]. VHL syndrome 
is characterized by cysts and tumors growing in specific parts of 
the organism [16-17]. It is considered a severe autosomal 
dominant genetic condition with inheritance of one person in 
over 35,000 [18]. The tumor injuries, which can be either benign 
or malign, are usually located in the retina, adrenal glands. 
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epididymis, central nervous system, kidneys and pancreas [19]. As 
a genetic disorder, VHL syndrome follows Knudson's two hit 
principle. A copy of the gene is mutated in the germ line, but the 
other gene copy still produces a functional protein. Complete 
protein inactivation appears during life due to somatic inactivation 
of the remaining functional copy [20] . On the contrary, mutations 
occurring during early fetal formation result in unsuccessful 
development [21]. The pVHL gene has 11,213 base pairs 
including three exons [18] and the final transcript is a protein 
commonly present in two isoforms: pVHLSO and pVHL19, of 213 
and 160 residues respectively. Neither isoform contains a known 
enzymatic domain, but rather appears to serve as a multipurpose 
adapter protein engaging in multiple protein-protein interactions 
[22]. pVHL structure is organized in an ot- and fi-domain and its 
stability was demonstrated to be ensured by direct interaction with 
other proteins such as Elongins B and C [23]. Both Elongin B and 
C are also required for the best characterized function of pVHL, 
the ubiquitination dependent degradation of Hypoxia Inducible 
Factor (HIF) via the proteasome [24]. However, pVHL is 
considered a multipurpose protein due to its high number of 
known interactors. At the time of writing, the IntAct database [2S] 
presents more than 200 diflFerent interaction partners, with some of 
them competing for the same Elongin binding site. Indeed, pVHL 
was found in different cellular compartments and seems to be 
involved in many difiTerent cellular processes such as apoptosis, cell 
proliferation, survival and motility [26]. Considering the huge 
number of interactors and multiple cellular localizations, many 
different functions have been described or hypothesized, such as 
regulation of cytoplasmic microtubules during mitosis [27] and 
endothelial extracellular matrix deposition [28]. On the other 
hand, considering the huge number of players involved in VHL 
syndrome and the lack of reliable kinetic data, a PN based 
approach may be a preferable option for an entire VHL pathway 
simulation. 

Petri Net for Interaction Pathways 

Since their invention, by Carl Adam Petri in the early sixties, 
PNs were mosdy used to describe technical systems, but later the 
utility in describing biological and biochemical functions has also 
been demonstrated [29]. PNs were successfully used in many 
studies to describe biological networks [30], such as the regulation 
and etiopathology in human Duchenne Muscular Dystrophy [3 1] 
and the hypoxia response network [12]. PNs are qualitative 
mathematical models that can graphically represent many object 
types, not only metabolites but also different protein states and are 
useful to simulate networks where not only metabolites are 
involved. Indeed, PNs can be a powerful tool to study all 
concurrent interactions in a specific pathway, even if the proteins 
or kinetics are not well-known. Due to the large number of 
different pVHL functions involved in VHL disease progression, we 
decided to extend the PN based analysis of [12] increasing the 
number of considered protein-protein interactions. We generated 
a novel manually curated PN model of the entire VHL regulation 
system collecting data from the literature and including the 
signaling pathways and glucidic metabolism. In order to build a 
realistic network, literature from both biochemical experiments 
and in silico predictions were used as source. It was decided to build 
a PN with only confirmed pVHL interactions whose function was 
also known. The resulting PN was validated using an analysis of 
specific properties as suggested by previous studies using the same 
method [29] . After validating the PN structure, in silico knock outs 
of specific proteins were done in order to observe the different 
network behaviors and the resulting biological effect. 



Methods 

The network was designed in the Snoopy PN framework 
(version 2, revision 1.13) [32], respecting the mathematical PN 
formalism as described in [12-33]. PN were demonstrated to be 
useful in describing discrete and concurrent processes in a simple 
graphical representation [30] and have been used to describe 
biomedical processes due to their capacity of representing 
sequential steps in a process. PN modeling methods are actively 
used to describe, simulate, analyze, and predict the behavior of 
biological systems. The Snoopy PN framework provides an 
extensible multi-platform framework to design, animate, and 
simulate Petri nets [32]. We chose Snoopy to facilitate future 
extensions of the VHL pathway presented here. Among different 
available PN types a standard PN was chosen to limit the number 
of variables. Both Charlie and PInA analyzers were used for PN 
analysis and validation (34). Further, in silico knock out experi- 
ments were used to test the biological reliability of the model. 
Structural model vaKdation was made by analysis of the T- 
invariants to demonstrate whether the system was covered by T- 
invariants and to confirm the biological meaning of each invariant. 
The use of T- and P-invariants is given by their own properties: 
they are a set (of transitions or places, respectively) that allow the 
reproduction of the same state after n transformations. A P- 
invariant represents a set of places where the number of tokens is 
constant and independent on the firing rate. A T-invariant instead 
represents a set of transitions that cyclically comes back to show 
the same initial set. Biologically a P invariant can represent the 
process of regulating a protein, whereas T invariants can represent 
cycUcal biochemical transformations such as metabolic reactions. 
To this end, the computed invariants were grouped in Maximal 
Common Transition Sets (MCTS) and Clusters, the former based 
on occurrence of specific sets of transition inside the various T- 
invariants, and the latter based on similarities between T- 
invariants. Diff(;r(;nt numbers of clusters will be defined depending 
on the resulting square matrix. Where MCTS create disjunctive 
nets, Clusters merge together similar T-invariants. Behavioral 
validation was made by selectively deleting tokens inside the 
model, imitating possible biological disruptions such as disease- 
causing mutations. The resulting network behavior was compared 
to what is reported in the literature. Total runtime for invariants 
computation were less than ten seconds on a mainstream Linux 
x86 workstation. Literature sources used to build the model are 
reported in Table SI. The Snoopy framework for PN construc- 
tion, Charlie and PInA tools for analysis are available at the 
website (URL: http://www-dssz.informatik.tu-cottbus.de/DSSZ/ 
Software). Finally, the model was used to simulate the network 
behavior through visual inspection of both token movement and 
accumulation in specific parts of the network. For a visual 
explanation of token movement in a PN refer to Video S 1 . 

Model Availability 

The resulting VHL disease PN model is available in File S 1 . 

Results 

Notations and Assumptions 

The PN built here focuses on pVHL interactions that were 
already proven by biochemical experiments and reported in the 
Uterature. We chose to model a realistic VHL disease pathways 
based on confirmed literature data, including all known VHL 
functions, VHL related signal pathway and glucidic metabolism. 
AU bibUographic sources used to design the model are presented in 
Table SI. The final PN is composed of 323 places and 238 
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transitions, connected by 801 arcs. Tables SI and S2 show all 
places and transitions and the related biological correspondence. 
Places are mainly proteins and enzymes, while some represent 
DNA or small molecular substrates such as glucose and cofactors 
(e.g. ATP). Notation for both pre- and post-places and their 
biological meaning are explained in Table SI. In a few cases, 
places are used to represent a whole group of changes generated 
by DNA transcription, (e.g. p_32 and p_33 or Et_efn and 
Et_efI2). Transitions instead symbolize complex formation be- 
tween two proteins or post-translational modifications. Output 
transitions stand for degradation or movement to other parts of the 
cell or organism to complete their functions (e.g. degrad_l and 
degrad_2) whereas input transitions show the generation of a 
substrate or protein. In order to simplify the design of such a large 
network, we decided to use macro nodes to group reactions 
representing complex molecular pathways such as signaling 
pathways or secondary signal cascades. The whole process is 
merged into a single node with a given name to allow visual 
inspection only in case of need. From the top level all transitions 
can still be found in a hierarchical lower layout level. Logic nodes 
were used for places participating in many reactions throughout 
the network such as ATP and ADP (7 logical copies each) or NAD 
and NADH (4 logical copies each). A total nesting depth of two 
was chosen to model macro nodes. Special arcs were not used 
while we chose to model the permanent presence of some objects 
using double arcs (e.g. for elob, eloc and places standing for 
enzymatic activity). In case of proteins which are actively 
degraded, it was preferred to create an input transition simulating 
constant production (or synthesis) and an output for consumption. 
This is the case for pkcz2, Jadel, pVHL and HIF-la. As can be 
seen from Figures 1 to 3, which represent the entire model, two 
major nodes can be immediately identified: pVHL and vcb, the 
complex made by pVHL and the two elongins. Another relevant 
part is the glucicUc metaboUsm, modeled due to its hypoxia 
induced regulation. It is represented in detail in Figure 2. 

HIF-la Transcription Activity 

The HIF-la transcription factor stimulates proliferation of 
endothelial cells to create new blood vessels during localized or 
broad hypoxia. In human, it is present as three different paralogs: 
HIF-la, HIF-2q( and HIF-3a. The sequence is quite conserved 
between the former two, whereas the latter is slightly shorter and 
seems to ha\e completely different functions compared to the 
other two [33-34]. Both HIF-la and -2a stimulate DNA 
transcription but the exact products of this activity are still poorly 
understood. In our model, only HIF-la in vivo activity was 
considered. It cannot be excluded that other biological effects 
depend on the second paralog. Indeed, both have a pro- 
angiogenetic function and are degraded by pVHL via proline- 
directed hydroxylation. HIF is a heterodimer of HIF- 1 a and HIF- 
IP, the latter being also termed Aryl hydrocarbon Receptor 
Nuclear Translocator (ARNT). We started from the transcription 
activity of HIF due to its regulation is the most studied pVHL 
function. Our model, as expected from literature data, shows that 
HIF-la enters the nucleus when not degraded by pVHL. It 
subsequently binds HIF- 1 P to form the HIF heterocomplex which 
interacts with DNA. Our model correctly simulates the increased 
affinity of HIF towards DNA. Transcription is enhanced by some 
co-factors binding both subunits of HIF and other proteins such as 
p300, Creb and cjun. This takes place in a specific DNA promoter 
sequence termed Hypoxia Response Element (HRE). Further- 
more, during transcription some pro-angiogenic factors are 
produced: Vascular Endothelial Growth Factor (VEGF), Endo- 



theUn (ET) and Erythropoietin (EPO). AU described pathways are 
in agreement with previous observations reported in [35] . 

Metabolic Processes 

HIF-la transcription activity includes some proteins which are 
dependent on oxygen but involved in other pathways (e.g. 
oxidative metabolism) or completely independent (e.g. metallo- 
proteinase MTIMMP). Further, HIF-la stimulates production of 

proteins imolved in the glucidic pathway. The final product of the 
metabolism is adenosine triphosphate (ATP), a molecular form of 
energy, composed by adenosine, an adenine ring connected to a 
ribose sugar, and three phosphate moieties. When a phosphate 
moiety is hydrolyzed it releases energy, used by cells for enzymatic 
reactions. The glucidic metabolism is composed of glycolysis, 
Krebs cycle, glycogen formation and respiratory chain with ATP 
synthesis. Glucose is absorbed in cells by enzymatic glucose 
transporters (GLUT), which carry the molecule to the location 
inside the cell where the metabolism takes place [36]. There are 
many isoforms of these transporters: GLUTl is present in all cells 
and in particular in erythrocytic membranes, neurons and glia 
[37]. GLUT2, located in both liver and pancreatic beta cells, is 
characterized by low affinity for glucose, hence it requires a higher 
glucose concentration to be activated [38]. Right after eating, 
glucose concentration increases, thereby quickly activating them. 
GLUT2 stimulates production of insulin, a hormone regulating 
the plasmatic glucose concentration. Glucose plasmatic concen- 
tration can also increase due to an opposite pathway, originating 
from liver glycogen being decomposed into glucose and reaching 
systemic circulation. GLUT3 is mostiy present in neurons, 
whereas GLUT4 is the insulin activated transporter located in 
myocytes, adipocytes and cardiomyocytes [36-39]. In our model, 
we chose to exclude GLUT3 due to its specific role in neuronal 
cells. Glycolysis occurs in the cytoplasm and during this process 
each glucose molecule is phosphorylated, consuming two mole- 
cules of ATP, then divided into two smaller molecules. Further 
modifications of these two molecules result in new ATP 
production. The molecule obtained at the end of glycolysis is 
pyruvate, which can be again modified through three different 
pathways. It can be decarboxylated and linked to Co-enzyme A to 
form acetyl-Co-enzyme A. It can then be carboxylated to obtain 
oxalacetate, or transformed through lactate dehydrogenase into 
lactic acid. Pyruvate can also be generated by other metaboUc 
pathways, like protein or fatty acid disruption and amino-acid 
modifications. Acetyl-CoA and oxalacetate are the molecules used 
in the following glucidic metabolism process, the Krebs cycle, 
taking place in the mitochondrial matrix. The Krebs cycle starts 
with acetyl-CoA and oxalacetate merging to create citric acid, 
which continues undergoing modifications until oxalacetate is 
formed again. During the process some co-enzymes are modified. 
Decarboxylation of pyruvate to form acetyl-CoA already trans- 
forms a NAD+ (Nicotinamide Adenine Dinucleotide) in NADH 
(reduced form), afterwards obtaining one more of ATP, GTP, 
FADH2 (Flavin Adenine Dinucleotide) and three more NADH per 
pyruvate molecule entering the Krebs cycle. The redox co- 
enzymes are considered electron transporters. During metabolic 
reactions they reduce themselves and get electrons (and protons) to 
oxidize the substrate of the enzymatic reaction. Electrons taken 
during the glucose metabolism are then used in the respiratory 
chain taking place in the internal mitochondrial membrane. The 
respiratory chain consists in transporting electrons through 
enzymes called cytochromes and others co-enzymes, characterized 
by the capability to receive and donate electrons. NADH (FADH2) 
is oxidized again by cytochromes going back to the form of NAD 
(or FAD). Electrons gained through oxidation are used to reduce 



PLOS ONE I www.plosone.org 



3 



June 2014 I Volume 9 | Issue 6 | e96986 



Petri Net Model of the VHL Tumor Suppressor 




Figure 1. Top level model. The colors of some tokens were arbitrarily chosen to give a clearer identification of the central nodes (ATP, Vcb and 
oxygen) or for nodes involved in more reactions such as GSKSp. The group of nodes in the bottom left is not disconnected from the central body of 
the network thanks to the presence of logic nodes for ATP synthesis (t_97). 
doi:10.1371/journal.pone.0096986.g001 



half a molecule of oxygen into water, releasing more energy. The 
FADH2 and NADH redox chain establishes a chemical potential 
causing the push of protons outside the internal membrane 
towards the inter-membrane space, which stays between the 
mitochondrial inner and outer membrane. This also causes a 
higher concentration of protons outside the inner membrane. The 
resulting gradient causes the tendency of protons to enter the cell. 
The final step is ATP-synthetase, formed by a channel that allows 
protons to enter, pushed by the gradient, allowing the enzyme to 
change conformation and make its reaction. This kinetic energy is 
converted into ATP. In our model, the glycolytic and Krebs cycles 
were described in detail, represented at the hierarchical second 
level by the coarse transition Glycolysis. The respiratory chain was 
instead merged into a single node (t_97). We chose to represent 
creation and consumption of ATP in order to show the effects of 
lower and higher oxygen concentration on the network. On the 
other hand, oxygen consumption for ATP synthesis during the 
respiratory chain creates a flow of oxygen in the model. Oxygen is 



not the only connection between glucidic metabolism and hypoxia. 
Indeed, HIF- 1 01 transcription activity enhances the transcription of 
many GLUT isoforms (such as 1, 3 and 9) and the pyruvate 
dehydrogenase kinase, which determines the pyruvate dehydro- 
genase (PyrDH) inactivation and consequent Acetyl-CoA forma- 
tion from pyruvate. Finally, Lactate dehydrogenase is also 
produced, to ensure an alternative compound, creating energy 
needed for cell survival [36-40]. 

pVHL-dependent Processes 

Some interactors can bind pVHL in regions interacting with 
Elongin C. These are HuR, Nur77, p53 andjadel. Nur77 has a 
complex function and its role in pVHL tumor suppressor activity is 
still not entirely clear. Nur77 can bind pVHL, inhibiting Elongin 
binding while allowing HIF- la binding. Its transcription is 
stimulated by HIF- la itself, and pVHL-HIF-la-Nur77 complex 
formation stabilizes the transcription activity of HIF- la by 
inhibiting the pVHL-dependent degradation [41]. Another 
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Figure 2. Lower hierarchical PN levels. Pathways from the top level are grouped in macro-nodes (functional subordinated layer), in particular 
glucidic metabolism and various VHL functions. 
doi:1 0.1 371/journal.pone.0096986.g002 
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Nur77 function is the stimulation of proopiomelanocortin 
(POMC) transcription, which is a precursor for adrenocorticotro- 
pic hormone (ACTH) formation. This hormone has an important 
stress response function, stimulating Cortisol production and other 
neurotransmitters from the adrenal glands, to enhance the 
organism reaction to danger and stress stimuli e.g. increase of 
gluconeogenesis and muscle mass. An excess of this hormone can 
cause desensitization of its receptors for feedback down-regulation 
and thus muscular weakness, tiredness, hyperglycemia and 
osteoporosis [42]. p53 can bind to pVHL avoiding the degradation 
of this tumor suppressor. Instead, it stimulates the apoptotic signal 
cascade via the p300 co-activator, which stimulates production of 
proteins enhancing the cell programmed death. If p53 cannot bind 
pVHL, two more mechanisms are described in the model. One is 
its modification and degradation by Mdm2 and the other is the 
pVHL-independent degradation of HIF-lot. Interaction with 
Mdm2 is needed in both cases [43-44]. Jadel is a short-lived 
protein whose main function is to stimulate the phosphorylation- 
dependent degradation of |3-catenin. This is a subunit of the 
cadherin protein complex acting as an intracellular signal 
transducer in the Wnt signaling pathway. It seems that P-catenin 
is able to stop cell division via a contact-dependent inhibition 
signal, whereas in Wnt signaling it is also involved in proliferative 
transcription. When Wnt is not present, P-catenin can be 
phosphorylated by Glycogen Synthetase Kinase, type SP (GSKSP) 
in complex with APC (Adenomatous Polyposis Coli) and Axin. P- 
catenin can interact with Jadel and be only successfully degraded 
after this interaction [45] . Related functions are represented in the 
macro node Jade l_pat. GSK3P seems to be a protein involved in 
many different pathways. GSK3P is involved in Glycogen 
Synthetase deactivation and can even phosphorylate pVHL and 
HIF-la. In the case of HIF-la, it generates a pVHL-independent 
degradation pathway, where phosphorylation allows ubiquitina- 
tion, whereas in the case of pVHL, it inhibits pVHL stabilization 
of microtubules [46]. 

Structural Model Analysis 

Based also on previous obser\'ations of Heiner et al., [47], in 
2008 Grunwald et al., demonstrated that PN can be used to 
describe large and complex metabolic pathways [31]. They 
postulated the following set of minimal rules that a PN should 
satisfy to be considered biologically reliable: (1) the network should 
be entirely connected, (2) the network should be covered by T- 
invariants, and (3) each T-invariant and P-invariant should have a 
biological meaning. The model described here was tested with 
respect to what previously done by Grunwald and co-workers [3 1] 
and resulted to be covered by T-invariants, connected, homoge- 
neous and each place has a pre-transition and a post-transition. 
Transitions without pre- or post-places were used to simulate the 
system interface to the surroundings. The network is alive, in other 
words, it continues to work forever, with all transitions contrib- 
uting to the net behavior forever, and no dead transitions. The 
MCTS and Cluster analysis were used due to the large number of 
T- and P-invariants included in the model. Both methods are used 
in PN theory to reduce the complexity connected with such a large 
network and to reduce the errors connected with manual 
investigation. From the 238 transitions present at the beginning 
in the model, 393 T-invariants were computed without consider- 
ing 10 trivial invariants. The latter consist in a pair of transitions 
that usually represent a forward and backward reaction, such as 
the active and inactive state of a protein. Trivial invariants could 
be erased to reduce the dimension of the network without 
disturbing the overall system when the interest is focused on the 
steady state behavior [12]. T-invariants were grouped into 44 



Clusters using the Tanimoto coefficient with similarity threshold of 
65%, as described in [31]. Only 11 of these 44 comprised more 
than one T-invariant. The three biggest Clusters are C9, 
composed of 144 T-invariants, C8 of 72 and Cll of 64 T- 
invariants. Separation into clusters allows easier analysis of 
networks pathways represented by each T-invariant, since they 
are grouped by similarity, specifically the common transitions by 
which they are composed. T-invariants named in the text are 
shown in Table S3, while T-invariants grouped in C8, C9, CIO, 
Cll are explained in Table S4 and described as follows. 

Cluster C8 

Cluster C8 groups all transitions included in HIF-la pathways, 
including transcription, signaling cascades, degradation via pVHL, 
p53 and GSK3P, and e\entually the Krebs cycle. For the EPO 
signaling pathway, two transitions (t 35 and t_36) are not included 
which cause Jak activation and consequent Stat5 activation to 
stimulate DNA transcription. Matrix stability regulation is also 
part of the cluster due to the destabilization induced by HIF-la 
transcription of metallo-proteinase (MMP), transitions from t_134 
to t_140. The largest T-invariant in C8 is Inv_280 (93 transitions) 
while the smallest is Inv_377 (81 transitions). The differences 
between T-invariants show the possibility of alternative pathways 
inside the model. For example, the VEGF dependent signal 
cascade can proceed in three different ways: t_13, t_14 and t_15, 
which lead to the pathways being merged in the coarse nodes 
Vegf_path3, Vegf_path2 and Vegf_pathl, respectively. The 
occurrence rate in C8 is 24 transitions for each path. The 
Endothelin, VEGF and Erythropoietin pathways are not in 
conflict and occurring together. Disaggregation of the matrix via 
MMPs is present in 18 T-invariants, whereas inhibition of these 
proteins, i.e. matrix stabilization, is present in the remaining 54 
transitions. Regarding the Krebs cycle, 47 T-invariants have t_91, 
of which only 24 reach t_92 and t_93, representing the last three 
steps of the cycle: succinate to fumarate, fumarate to malatc, and 
malate to oxalacetate. AU the malate being produced is used to 
regenerate oxalacetate. Degradation of HIF-la occurs in any T- 
invariant of the cluster. The pVHL-dependent degradation of 
HIF-la is always present (transitions t_116 to t_119). In 19 T- 
invariants degradation takes place via p53 (t_191 to t_193) or, 
alternatively, via phosphorylation by GSK3P in another 17 T- 
invariants. Two of the three pathways can be present in the same 
T-invariant, as in Inv_227, where degradation via pVHL and 
degradation via p53 are both present. This was considered as the 
HIF- 1 a dependence on the lack of degradation by these proteins. 
AU three degradation pathways never appear in the same T- 
invariant. The p53 and GSK3P paths are never present together 
but each of them is accompanied by pVHL-dependent proteaso- 
mal degradation. Inv_377 lacks the EPO signaling pathway but is 
the only one in this cluster to have t_34, t_33 and t_37. These 
invariants have all input and output transitions. For example, 
t_202 the second input for pVHL, is present in only 18 invariants. 
Other inputs are t_98, always present, leading to formation of 
HIF-la and pVHL, t_192, producing p53 and t_216, representing 
other pyruvate generating metabolic pathways. The latter is also 
present in each invariant allowing formation of the pyruvate 
needed for Krebs cycle progression. 

Cluster C9 

Cluster C9 is the largest cluster in our model and includes 144 
T-invariants. It is characterized by complete EPO pathway 
abrogation which goes through formation of the Shc-Grb-Sos 
complex and the consequent mapk-dependent phosphorylation 
cascade. Transition t_127, representing EPO effects on oxygen 
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Table 1. List of Trivial T-invariants excluded from calculation with their associated biological meaning. 



Trivial T-invariants 


ID transitions 


Biological Meaning 


Tlnv_l 


t_99, t_100 


Glycongen Synthase regulation 


Tlnv_2 


t_101, t_102 


Pkb regulation 


Tlnv_3 


t_103, t_132 


GSK3P active-inactive state 


Tlnv_4 


t_l 74, t_222 


Par6 inactivatlon via aPKCt^2 


Tlnv_5 


t_l 77, t_208 


VHL binding to Spl 


Tlnv_6 


t_167, t_199 


Spl phosphorylation and dephosphorylation 


Tlnv_7 


t_0, t_2 


Hif transport in and out of nucleus 


Tlnv_8 


t_0, t_234 


Hlf inhibition via FIH 


Tlnv_9 


t_181, t_207 


Hur inhibition via VHL 


Tlnv_10 


t_231, t_232 


IGFR mRNA production and destruction 


doi:l 0.1 371 /journal.pone.0096986.t001 



production, is absent. In its place, t_35 and t_36 are considered, 
wliich are present in 72 T-invariants. In cluster C9, the largest T- 
invariants are Inv_278 and Inv_279 (74 transitions) while the 
shortest ones are Inv_101, Inv_105, Inv_144 and Inv_148 with 65 
transitions each. 

Cluster CIO 

Cluster CIO, composed of 52 T-invariants, is characterized by 
the presence of glycolysis between many transitions grouped in the 
cluster. This is also the cluster containing the most populated T- 
invariant of all computed 393 non-trivial T-invariants. This is 
Inv_245, including 101 transitions and covering almost half of the 
whole model. Cluster CIO also includes Inv_125, the shortest 
invariant of this model, composed by 85 transitions due to lack of 
the Krebs cycle. Another difference with the other three major 
clusters is that here both EPO paths are present, specifically, the 
Jak pathway belongs to 4 T-invariants and Shc-Grb-Sos is 
observed throughout the cluster. Vegf_pathl seems to be more 
common in this cluster, being present in 36 T-invariants, whereas 
the other two are present 12 times each. This time they are present 
even in the same invariant, as for Inv_60, Inv_129, Inv_172 and 
Inv_215, with both t_13 and t_15, and Inv_142, Inv_185 and 
Inv_228 with t_14 and t_15 and all subsequent signaling 
appearing at the same time. Despite glycolysis being present in 
all cluster invariants, the Krebs cycle appears only in 1 1 cases. 
p53-dependent degradation of HIF- lot occurs in 1 1 cases while the 
phosphorylation-dependent one appears in 13. An input transition 
has been added with respect to the other major clusters so far 
analyzed (i.e. t_69_eating) without which glycolysis could never 
take place. 

Cluster C1 1 

Cluster CI 1 is composed of 64 T-invariants. Only part of the 
EPO pathway is described here, with the major difference that the 
Krebs cycle is completely abrogated while Prolyl Hydroxylase type 
2 (PHD2) regulation by oxalacetate is included. HIF- lot interac- 
tion with Nur77 and transcription of VEGF by Spl are also 
present. t_80 (transformation of pyruvate in oxalacetate) is not 
present in the first 42 cluster T-invariants. Nur77 interaction with 
HIF- la is present only in 8 T-invariants, specifically Inv_87 to 
Inv_94. Spl transcription activity is appearing in twice the 
amount, including the same 8 invariants just mentioned. VEGF 
transcription via Spl activity is aPKC^2 phosphorylation depen- 
dent, which does however not appear in the cluster. When VEGF 



is synthesized, it is subsequently stabilized by Hur, followed by 
t_178 and Hur is recreated to allow other functions. Indeed, it is 
one of the few places without input transition but with a token that 
goes forward and backward again. Compared to the other clusters, 
C 1 1 also shows one less transition in the coarse PHD regulation 
node, specifically t_81, which shows the transformation of 
pyruvate by pyrDH into acetyl-Coenzyme A, needed for the 
Krebs cycle. The four clusters C-8 to C 12 are very similar to each 
other, as can be seen from the distance tree in Figure 4. They all 
contain the HIF- lot transcription activity and signaling pathways 
caused by EPO, VEGF and the HIF- la degradation options. 
They include the effects of other transcription activity products, 
like metallo-proteinase and pyruvate dehydrogenase kinase, which 
regulate the activation state of PyrDH. All include part of the 
glucidic metabolism but not Glycogen formation itself Other five 
clusters from C 1 2 to C 1 6 have a smaller number of T-invariants 
and fewer transitions present in each invariant. They do not 
include transcription activity but are only formed by the VEGF 
and glycolytic pathways. The information contents of these clusters 
turned out to be uninformative and their analysis was not 
included. The same applies to clusters composed by 1-3 T- 
invariants. Finally, some transitions are not present in the clusters 
and not listed in the T-invariants because trivial invariants were 
excluded from cluster analysis. These transitions are shown in 
Table 1 with their respective biological meaning. 

MCTS Analysis 

Another way to group invariants is by the amount of single 
transitions present in them. Maximal common transition set 
(MCTS) analysis provides a PN decomposition into non-overlap- 
ping subnets, sharing parts of the same T-invariants [29] . In a 
biochemical network, MCTS could be interpreted as enzyme 
subsets operating together under steady state conditions, comput- 
ed based on the support of a T-invariant. MCTS computation 
does not consider stoichiometric relations, describing exclusively 
sets of reactions present in a maximal number of T-invariants 
resulting shared by different signaling pathways [48] . A total of 40 
non-trivial MCTS were identified, with results and related 
biological means shown in Table 2. Some transitions do not 
belong to a non-trivial MCTS, because their occurrence has no 
similarity with other transitions and they create separate MCTS 
(specificaUy: t_69, t_82, t_91, t_94, t_98, t_114, t_116, t_120, 
t_121, t_122, t_179, t_202, t_209, t_212, t_216, t_225 andt_229). 
MCTS define transitions that always take place together, but are 
not necessarily connected, thus representing disjunct building 



PLOS ONE I www.plosone.org 



8 



June 2014 | Volume 9 | Issue 6 | e96986 



Petri Net Model of the VHL Tumor Suppressor 



44 43 42 41 40 39 38 37 36 35 34 33 32 31 30 29 28 27 26 25 24 23 22 21 2019 18 17 16 15 14 13 12 1 1 10 9 8 



lj 



7 6 5 4 3 2 1 



Figure 4. PinA Distance Matrix clustering, using Tanimoto coefficient and 65% tKireshold of. The numbers indicates clusters. In C8, 
C9, CIO, C11 are highlighted a red square. 
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blocks constituting the network. Considering both analyses, a table 
was automatically built in PInA [29] showing a correlation 
between clusters and MCTS. Transitions (t) or MCTS (M) are 
compared to evaluate how many T-invariants clusters cover the 
selected M or t (if the transition is not already part of the MCTS, 
as listed above). The more covered a transition or set is, the more 
central it could be considered for the network behavior. Recently, 
a network coarsening method based on abstract dependent 
transition sets (ADT) was presented [49] . It is formulated without 
the requirement of pre-computation of the T-invariants and is a 
tool commonly used for the decomposition of large biochemical 
networks into smaller subnets. Due to the manually designed 
nature of our model, we preferred to maintain a logic hierarchy 
based on metabolic pathways in order to maintain the network 
centered on pVHL and its interaction. The MCTS calculation 
results shows that the most covered set by cluster T-invariants is 
M20 with 358 T-invariants covering all transitions in the set, 
indicating that this MCTS corresponds to more T-invariants than 
the others. AH transition sets are an important link to the others, as 
tokens pass through these transitions more often. A transition not 
present in any set but most covered by T-invariants is t_98, which 
is also the most frequently occurring transition, see Figure 5. The 
10 most occurring transitions are listed in the Table 3. 

P-invariant Analysis 

Although the network is not covered by P-invariants, it has 130 
P-invariants. 47 of these are trivial P-invariants, comprising a 
single place, connected with double arcs to imitate an activator arc 
function. Another object represented with double arcs is the 
enzymatic activity catalyzing a reaction and immediately going 
back to the steady state. P-invariants show places or sets of places 
where token numbers always remain equal and do not move 
outside the subnetwork induced by the P-invariant in the initial 
marking. In other words, they do not grow nor diminish. The 
remaining P-invariants are mostly located in signal transduction 
pathways, such as situations in which a protein is sequestered from 
its function and then goes back after a second reactivation 
mechanism. This scenario is present in p_41, p_42 and p_45 
located in invariant P_58. It is important to notice that ATP and 
ADP, as well as NAD and NADH, are modeled as P-invariants. 
P_90, P_91 and t_97 are able to transform ATP and ADP. More 
in general, all energy consuming transitions are considered to be 
backward transitions of invariants. Invariants not related to signal 
transduction are places located in the Hur system, where Hur is 
removed from its function by pVHL. This is a good approxima- 
tion for sequential modifications that momentarily activate 
proteins. Afterwards, Hur can go back and stabilize VEGF to 
increase its transcription activity. 



In Silico Knock Out Experiments 

The previously described clustering and MCTS analysis for T- 
invariants allowed us to identify the most common transitions and 
to understand which transitions can be depleted in our knock out 
experiments in order to get the most important biological effect. 
The knock out experiments were performed erasing selected 
transitions or tokens and observing which transitions or MCTS 
become inactivated. Considering our results and the literature, we 
decided to knock out the following pathway elements: (i) pVHL, (ii) 
HIFlot alone and with Spl, (in) t_98, (iv) PHD2, (v) MCTSl, (vi) 
t_97 and (vii) GSKSp. In the following, we describe the effect of 
each knock out scenario on our model. 

(i) pVHL knock out. Degradation of HIF-la is not 
completely depleted due to presence of both p53- and GSK3|3- 
dependent alternative degradation pathways. AU other processes 
usually inhibited by pVHL take place in an uncontrolled way, 
including creation of VEGF via Spl transcription activity and 
increased matrix regulation due to lack of fibronectin crosslinking. 
Hur resulted constantly activated and nur77 can stimulate 
synthesis of Proopiomelanocortin, precursor for the Adrenocorti- 
cotropic hormone. Card9 increases release of tumor necrosis 
factor, and NF-kB when not inhibited by p VHL. Instead, Jade 1 is 
unable to survive long enough to inhibit (5— catenin, generating a 
proliferation signal with Wnt. Lactic acid is also not produced due 
to LDH enzyme production being HIF-la transcription activity 
dependent. 

(ii) HIF-la knock out. VEGF is still created thanks to Spl, 
thus oxygen is still generated even if in lower proportion. If HIF- 
la and Spl are both knocked out at the same time, oxygen is 
quickly consumed and the metabolism is soon unable to proceed. 
Lactic acid is not produced due to LDH enzyme production being 
HIF-la transcription activity dependent. Glycolysis and glycogen 
are produced normally and the metabolism is not inhibited by 
PyrDH negative regulation and lactic acid formation. Since pVHL 
is present, other tumor suppressor activities are enabled, except for 
proteasomal degradation of HIF-la due to the substrate being 
non-existent. 

(iii) HIF-la and pVHL double knock out. This generates a 
situation where the metabolism is normal but oxygen regeneration 
is less productive, with only Spl acting for transcription. Due to 
absence of pVHL, all proliferation-stimulating processes are 
active, causing an unbalanced consumption of resources. Our 
model shows that this condition is compatible with cell growth and 
multiplication, but new blood vessel generation is consistently 
slower and glucidic metabolism appears principally based on the 
glycolysis reaction. Similar activity reduction applies to both tight 
junction and cellular external matrix (ECM) pathway regulation. It 
cannot be excluded that some observed effects could be mitigated 
by both HIF-2a and HIF-3a activity in vivo. 
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Table 2. List of IVICTS and transitions from PInA. 



MCTS ID Transitions 



MCTS 1 (Ml) 


t_a t_190, t_191, t_192; 




MCTS 2 (M2) 


t 1, t 3, t 4, t 5, t 6, t 7, t 8, t 9, t 10, t 11, t 12, t 32, t 51, t 52, t 53, t 54, t 55, 
t 65, t 66, t 67, t 79, t 83, t 123, t 129, t 138, t 215; 


t 56, t 57, t 58, t 59, t 60, t 61, t 62, t 63, t 64, 


MCTS 3 (M3) 


t 2, t 99, t 100, t 101, t 102, t 103, t 132, t 167, t 174, t 181, t 199, 177, t 207, t 


208, t 222, t 232, t 234; 


MCTS 4 (M4) 


t 1 3, t 1 6, t 1 7, t 1 8, t 1 9, t 20, t 1 24, t 1 28; 




MCTS 5 (MS) 


t 14, t 21, t 23, t 24, t 25, t 26, t 27, t 31, t 126; 




MCTS 6 (M6) 


t 15, t 22, t 28, t 29, t 30, t 125, t 142; 




MCTS 7 (M7) 


t_33, t_34; 




MCTS 8 (M8) 


t_35, t_36; 




MCTS 9 (M9) 


t_37, t_l 30; 




MCTS 10 (MIO) 


t 38, t 39, t 40, t 41, t 42, t 43, t 44, t 45, t 46, t 47, t 48, t 49, t 50, t 127; 




MCTS 11 (Mil) 


t 68, t 70, t 71, t 72, t 73, t 74, t 75, t 76, t 77, t 78, t 237; 




MCTS 12 (Ml 2) 


t_80, t_l 1 1 , t_1 1 2; 




MCTS 13 (Ml 3) 


t 81, t 84, t 85, t 86, t 87, t 88, t 89, t 90, t 131; 




MCTS 14 (Ml 4) 


t_92, t_93; 




MCTS 15 (Ml 5) 


t_95, t_104, t_141; 




MCTS 16 (Ml 6) 


t_96, t_105, t_218; 




MCTS 17 (Ml 7) 


t_97, t_224; 




MCTS 18 (Ml 8) 


t_106, t_107, t_108; 




MCTS 19 (Ml 9) 


t_109, t_110, t_133; 




MCTS 20 (M20) 


t_l 13, t_l 15, t_1 17, t_l 18, t_1 19; 




MCTS 21 (M21) 


t_134, t_135, t_136, t_137; 




MCTS 22 (M22) 


t_139, t_140, t_217; 




MCTS 23 (M23) 


t_l 43, t_227, t_228; 




MCTS 24 (M24) 


t_144, t_145, t_146; 




MCTS 25 (M25) 


t_147, t_148, t_149; 




MCTS 26 (M26) 


t_150, t_151, t_219; 




MCTS 27 (M27) 


t_l 52, t_l 53, t_l 54, t_220; 




MCTS 28 (M28) 


t_155, t_156, t_157, t_158, t_159, t_160, t_221, t_226; 




MCTS 29 (M29) 


t_161, t_163, t_164, t_166, t_213; 




MCTS 30 (M30) 


t_162, t_165, t_214; 




MCTS 31 (M31) 


t_168, t_169, t_201, t_236; 




Iv\K- 1 J DZ \IV\jZ) 


t 1 7n t 171 t 177 t 17^ t ll^- 
l_ 1 /U, I_ 1 / 1 , I 1 / Z, I_ 1 / D, t ZZj, 




MCTS 33 (M33) 


t_175, t_176, t_180, t_230; 




MCTS 34 (M34) 


t_178, t_203; 




MCTS 35 (M35) 


t_182, t_183, t_231, t_233; 




MCTS 36 (M36) 


t_184, t_185, t_186, t_187, t_188, t_189, t_198; 




MCTS 37 (M37) 


t_193, t_194; 




MCTS 38 (M38) 


t_195, t_196, t_197, t_200; 




MCTS 39 (M39) 


t_204, t_205, t_206, t_235; 




MCTS 40 (M40) 


t_210, t_211; 





doi:1 0.1 371 /journal.pone.0096986.t002 



(iv) PHD2 knock out. The protein is involved in pVHL 
mediated and oxygen dependent degradation of HIF- 1 ot. Further, 
PHD2 is involved in hydroxylation of tlie RNA polymerase II 
subunit Rpbl to allow its translocation to less chromatin- 
concentrated areas of the nucleus. When it is knocked out, HIF- 
la degradation can continue via alternative pathways as seen in 
the pVHL knock out experiment and there is more RNA 



polymerase II activity, even if rpb7 can stiU be inactivated by 
pVHL. 

(v) MCTSl knock out. MCTSl groups some reactions 
involved in the HIF- lot p53-dependent degradation pathway 
(Table 2). To perform this knock out, we erased the necessary 
token in mdm2, making the precondition insufficient to enable the 
MCTS transitions. p53 is not degraded and can continue its 
proapoptotic signal. On the other hand, a HIF- la degradation 
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Table 3. Ranking of the 10 most occurring transitions with biological meaning and percentage of occurrence. 



Rank 


Transitions 


Biologicsl meBninci 


Occurrence % 


1 


t_98 


Input transition for Hif and VHL 


95.165 


2 


t_116 


Interaction of VHL with Elongin B and C 


94.148 


3 


t_113 


Activation by oxygen of ARD 


94.094 


4 


t_115 


Acetylatlon and hydroxilatlon of Hif 


94.094 


5 


t_117 


Interaction of complex Vcb with Cu2 


94.094 


6 


t_118 


Interaction of complex Vcb with modified Hif 


94.094 


7 


t_119 


Degradation VHL dependent of Hif 


94.094 


8 


t_97 


ATP formation 


89.059 


9 


t_224 


Water Output transition 


89.059 


10 


t_82 


Pyruvate Dehydrogenase inactivation 


88.041 



doi:l 0.1 371 /journal.pone.0096986.t003 



mechanism is also knocked out resulting in an increased HIF- la 
transcription activity. 

(vi) t_97 knock out. This is the ATPase transition, allowing 
the model to imitate oxygen consumption for ATP synthesis. If this 
transition is inactive, oxygen accumulates infinitely and ATP is not 
regenerated after few simulation steps. At the beginning, ATP is 
formed during the first step of glycolysis but afterwards it is 
consumed again. At some point, these reactions do not have any 
ATP available to allow the system to re-balance the consumed 
ATP. After few simulation steps, oxygen reaches a high level due 
to slower consumption in the PHD2 regulation process. Biolog- 
ically, this means that the metabolism stops and the cell is not able 
to create energy to survive. There is no accumulation other than 
glucose in the model. A few oxygen creation processes are blocked 
as well due to absence of ATP, e.g. t_15, t_41 and t_57. 

(vii) GSK3P knock out. This enzyme is involved in negative 
glycogen synthetase (GS) regulation and is inactivated when 
phosphorylated. When GSK3P is knocked out, glycogen is 
continuously produced due to the enzyme remaining in an active 
state. In a real organism there are alternative forms of GSK3P 
which can inactivate GS, hence the effect will be less sharp. 
GSKSP is also involved in the degradation of HIF- lot, causing its 
phosphorylation and following ubiquitination. It is also involved in 
the degradation of P-catenin, where it is responsible for primary 
phosphorylation. If knocked out, even if Jade 1 can be stabilized by 
pVHL, the effect will be similar to a knock out of Jade 1, where P- 
catenin is free to continue proliferation stimulating transcription 
activity. 



Discussion 

We started from a core model of hypoxia response [12] and 
extended the original network with functional data derived from 
the literature in order to represent a complete description of the 
pVHL interaction pathway according to current knowledge. VHL 
syndrome is characterized by the formation of tumors and cysts 
affecting different organism districts and tissues. Indeed, pVHL is 
a tumor suppressor whose functions are connected to inhibition of 
proliferation and survival, growth and stability of extracellular 
matrix and microtubules, as well as cell polarity and migration. 
The IntAct database reports more than 200 suspected pVHL 
interactors and for most of them interaction and function details 
remain largely unknown. We chose to model the pVHL 
interactions in a credible cellular context with many protein 
activities occurring at the same time. The main idea was to create 
a novel manually curated PN description of the entire VHL 
disease pathway, including glucidic metabolism and signaling 
pathways. The model was designed as a standard PN and is 
composed of 238 transitions and 323 places, connected by 801 
edges. A biologically realistic PN model needs to be covered by T- 
invariants, meaning each transition in the model has to be 
included in a T-invariant, and each invariant needs to have a 
biological meaning [31-47]. We used the T-invariant analysis to 
validate the reliability of the model. We computed a total of 393 
T-invariants, plus 10 trivial invariants, which were excluded from 
analysis. These were grouped into 44 Clusters and, through use of 
T-invariants, transitions were grouped into 40 MOTS. The model 
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Figure 5. Transitions occurrence T-invariants. Transitions are ordered by name and t_98 Is highlighted In red. 
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obtained is connert(;d, covered by T-invariants with each invariant 
holding a biological meaning. MOTS analysis was used to identify 
the most frequent crucial transitions occurring in the model. This 
specific subset was further used to plan in silico knock out 
experiments and for the model validation and analysis of expected 
biological behavior. The model was then used to perform in silico 
knock out experiments inactivating specific transitions during 
qualitative network analysis. Our results showed that the model is 
able to represent important transitions reflecting real biological 
outcomes, i.e. transitions involving species such as oxygen or ATP 
are correctly inactivated under cc'rtain circumstancx's as cxpc'cted 
from the bibliographic data. Biological energ^'-related reac:tions 
(e.g. ATP production from ADP) were modeled as P-invariants. 
Although the network is intentionally not covered by P-invariants, 
P-invariant analysis was used to verify all modeled energy 
consuming transitions. Both the ATP and NADH balances 
appeared constant during the simulation, with irrelevant P- 
invariants located in the Hur system. This approximation was 
used to verify the Hur-dependent regulation of VEGF, with results 
in accordance with [50]. The specific pVHL knock out suggests 
that this protein alone is not sufficient for complete HIF-la 
inactivation. Indeed, other concurrent HIF-la degradation 
pathways promote a sort of cell cycle regulation backup. On the 
contrary, simple deletion of pVHL turned out to be sufficient to 
increase all its other inhibitory functions, showing similar effc-cts to 
pathological VHL symptoms. Indeed, ECM destabilization 
increases cell migration to other areas, promoting metastasis 
outbreak in case of tumor cells. Further, pVHL-dependent 
inhibition of tight junction formation by aPKCpi participates in 
an easier cellular detachment. The interactions of Nur7 7 could be 
considered a good example for pathological effects. It is a 
stimulator of Proopiomelanocortin production, a precursor for 
the Adrenocorticotropic hormone. If excessively released, it 
promotes an overproduction of adrenergic neurotransmitters by 
adrenal glands. Coming at clinical condition known as Gushing 
syndrome. On the very long term, Nur77 deregulation is known to 
cause tumors of the pituitary and adrenal glands [51], [52]. This 
happens in pheochromocytoma, which is one of the main VHL 
disease manifestations. We speculate that continuous VEGF 
transcription, even in situations where HIF-la (but not Spl) is 
knocked out, could be the explanation for clinical studies where 
VEGF-targeting drugs have turned out to be effective in kidney 
cancer treatment as reported in [53]. Although we used only 
confirmed data from the hterature, Nur77 may be involved in 
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